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1. Motivation and Objective 

The first objective is to evaluate current two-equation and second order closure 
turbulence models using available direct numerical simulations and experiments, 
and to identify the models which represent the state of the art in turbulence mod- 

eling. 

The second objective is to study the near-wall behavior of turbulence, and to 
develop reliable models for an engineering calculation of turbulence and transition. 
The third objective is to develop a two-scale model for compressible turbulence. 

2. Work Accomplished 

2.1 Evaluation of two-equation models (N.J. Lang and T.-H. Shih 1 ) 
Twelve two-equation models including k - e, k — r , k- to and q - to models, 
have been evaluated using a common flow solver code (GENMIX) for wall bounded 
turbulent flows. For each model, calculations were carried out for two-dimensional, 
fully developed channel and flat plate boundary layer flows. These flows are ap- 
pealing for model testing because they are simple and self-similar, yet demonstrate 
important features of wall bounded turbulent shear flows. In addition, we can 
compare these calculations with Direct Numerical Simulations (DNS). A list of the 
models which were tested are shown in the table below: 


Ch 

Chien 2 1982 

fc — e 

Sh 

Shih 3 1991 

k - e 

LB 

Lam and Bremhorst 4 1981 

k - e 

NH i 

Nagano and Hishida 5 1987 

k — e 

NT 

Nagano and Tagawa 6 1990 

k — 6 

LS 

Launder and Sharma 7 1974 

k - € 

JL 

Jones and Launder^ 1973 

k — e 

MK 

Myong and Kasagi 9 1988 

k - e 

YS 

Yang and Shih lu 1991 

k — e 

WI1 

Wilcox 11 1984 

k — to 

WI2 

Wilcox 12 1991 

k — to 

SAA 

Speziale, Abid and Anderson 13 1990 

k - t 

Co 

Coakley 14 1983 

q- u) 


Two dimensional channel flow calculations were made at Re r = 180 and Re r = 
395. These cases were compared with DNS data of Kim et al 15 . Calculations for 
the two-dimensional flat plate boundary layer flow at Re$ = 1410 were compared 
with DNS data of Spalart 16 . Some flat plate boundary layer comparisons were made 
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between the experimental data of Klebanoff 17 at Re$ = 7700 and solutions of the 
various models. The detailed results are reported in NASA TM 105237, 1991. 

An important criterion for two-equation model comparisons is not only how well 
the model predicts mean velocity and shear stress, but also the turbulent kinetic en- 
ergy and dissipation (or specific dissipation) rate. These predictions should provide 
appropriate turbulent velocity and length scales so that the model can be applied 
to more complex flows for which a simple mixing length model often fails. Some 
researchers maintain that it is not critical whether or not the turbulent kinetic en- 
ergy and the turbulent length scale are predicted with great accuracy. However, 
one may imagine that a two-equation model making unreasonable turbulent veloc- 
ity and length scale predictions would be very questionable when applied to more 
general flows. A model which accurately predicts the shear stress and mean veloc- 
ity does not imply that it has correctly modeled the turbulent kinetic energy and 
length scale. In fact, if the turbulent kinetic energy is incorrect, then the length 
scale must also be incorrect to compensate for the error in the turbulent kinetic 
energy. For this case, two wrongs are making a right. This warrants some caution 
when computing flows for other geometries. 

Comments on two-equation models: 

It is clear that the JL, LS, WI1 and WI2 models underpredict the near wall 
turbulent kinetic energy compared to the other models. 

The standard k - e model has been proven to provide good results in the high 
Reynolds number range. It is therefore attractive for a near wall k - e turbulence 
model to approach the standard k — e model away from the wall. The LB, LS and 
VS models are the only k — e models in this study which possess this characteristic. 

Because boundary layer and channel flows are self-similar, the solutions should 
be independent of the initial conditions. However, some of the models (SAA, Co, 
and LB) have difficulty when the initial conditions contain large gradients. The Co 
Model is particularly dependent on the initial conditions. Even slight perturbations 
of the initial conditions will yield noticeably different solutions with this model. 

JL, LS, WI1 and WI2 are the only models which do not contain y + . Damping 
functions which contain y + are not desirable because y + is erroneous near flow 
separations and not well defined near complicated geometries. Unfortunately, these 
are the same models which poorly predict the near wall turbulent quantities. 

The Wilcox models (WIl and WI2) suffer from a numerically awkward boundary 
condition for u at the wall: 

w -> W 88 y+ 0 

We cannot define u; at y + = 0. We have tried two ways to approximate u as y + 
approaches 0. First, we chose a large number for u wa ii and, second, we used an 
asymptotic u) wa ii = Test cases showed that the solution does not converge 

as u wa u — ♦ oo when y + — ♦ 0 for either model. In addition, both Wilcox models 
underpredict the turbulent kinetic energy peak value for both boundary layer and 
channel flows. 
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model the new unknowns Ty, II y and e y . At the level of the second order closure, 
these new unknowns are usually modeled with algebraic equations in terms of the 
second moments and the mean quantities (with the exception of the trace e k k 

which is modeled with a transport equation). 

In a general turbulent shear flow with moderate inhomogeneity, the turbulent 
diffusion terms in the second moment equations are usually smaller than the ot er 
terms However, the pressure-strain rate and dissipation rate tensors are always 
among the leading terms. Therefore, the performance of modeled equations largely 
depends on the models of pressure-strain rate tensor and dissipation rate tensor. 

In this study, we only concentrate on the models of the pressure-strain rate tensor 
and the dissipation rate tensor for the velocity field: Hy - ey, which are modeled 


as 




where 11^ are respectively referred to as the rapid term and the return-to- 
isotropy term. 

Models for the rapid term n§* 

Launder, Reece and Ro di (LRR ): 18 


n 


+ 


where C 2 — 0.4, and 


=- = 0 . 2 Sij 4 - (bikSjk + bjkSik - -zbijbkiSki) 

2 q 2 22 6 

2 {bikftjk + bjk^ik) 


( 2 . 2 . 1 ) 


22 


_ UiUj 1 c. . 
bij - -== 


Sij-jiUij + Uj'i), 

n« = 

This model is linear in the Reynolds-stress. It contains only one model constant C 2 . 
This model satisfies the conventional model constraints 34 . It is the most general 
form at the level of Unear dependence on the Reynolds-stress. However, as Lumley 
pointed out, this model may violate readability as the turbulence approaches a 

two-component state. 

Spasiale. Sarkar and G atski (SSG ); 19 


2 q 1 4 2 q 2 


+ bikSjk + b jk Sik - | SijbkiSki ) 

4 ~r(bik^jk 4 " bjk^ik) 

4 


( 2 . 2 . 2 ) 
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where, 

CV = p = -whUij 

C{ = 1.8 C 3 * = 1.3, C 4 = 1.25, C 5 = 0.4 

This model is quasi-linear in the Reynolds-stress, because the coefficients in the 
first two terms are not constant, but depend on the invariant of the Reynolds- 
stress tensor and the production P, This model contains four model constants 
(Cf, £3,04,(75), therefore one may imagine that it will be difficult to correctly 
calibrate them. In addition, this model does not satisfy the normalization constraint 
which is one of the basic model constraints 34 . If we impose this constraint, then 
the four coefficients will reduce to only one, and this model will reduce to the LRR 
model. Finally, like the LRR model, the SSG model may also violate realizability. 
Fu. Launder and Tselepidakis (FLT): 20 

n^ } 2 

— =- = 0.2 Sij + 0.3(bikSjk +bjkSik — - SijbkiSki ) 

2 q 2 J 

1.3 

H {bik^ljk "I" bjk^tik') 
o 

+ 0 . 2 ( 6 ?, Sji + 6 ?, Si, - 2b kj buSki - 36 iy 6*,S fc , ) ( 22 ' 3) 

+ 0.2(6 + 6?,fij,) 

+ r[4bl n (b ik Cljk + bjkttik) 

"I" 126 rn ,6 n j(6 T7 ifcflnfc d" 6f,fcflmfc)] 
where r = 0.7, 6?^ = 6^6*^. 

This model is cubic in the Reynolds-stress. The final form selected contains 
one model constant. This model only satisfies a part of the realizability condition, 
corresponding to a two-component state of turbulence. However, when a scalar 
field is involved, this model cannot satisfy Schwarz’ inequality between velocity and 
temperature. This part of realizability is sometimes called joint realizability. 

Shih and Lumlev (SLli 21 

n< 1} 2 

===- — 0.2 S{j ”f" 3oi5 {bikSjk d- bjkSik ~ bijbkiSki ) 

2 q 2 A 

+ ^(2-7as)(bikCljk + bjkftik) ( 2 - 2 - 4 ) 

+ 0.2(6?, Sj, + &?, Su - 2bkjbuSki - SbijbuSki) 

+ 0.2(6?,flj, + tijiQ.ii) 

where, 

ctj = — (1 0-8F 1 / 2 ), F = Id- 9bijbjkbki — ^ bijbij 

This model is quasi-quadratic in the Reynolds-stress, because the model coefficient 
<15 is a function of the invariants of the Reynolds-stress tensor. We emphasize that 
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this model is obtained from a more general form of the expressioh than the FLT, 
and satisfies both the two component condition and Schwarz’ inequality between 
the velocity and scalar fields. In addition, the final form is simpler than the model 

of FLT. 

Shih and Mansour (SM): 22 

nW 2 

-±Lr = 0.2Sij + Za b (bikS jk + bjkSik - -zS^buSu) 

2q 2 6 

+ \(2- 7 a 5 )(b ik n jk + b jk n ik ) (2-2.5) 

o 

+ 0.2(6?,5ji + b^Su - 2b kj buS k i - 3bijb kl S k i) 

+ 0.2(6|jfljj + tijidii) 


where, a 5 = ^{1 + 3.5[1 - (1 - •F’) 1 ^]}- _ 

This model has the same form as the SL model. It was derived in a different way 
and contains a different model coefficient ct 5 which was calibrated from one of the 
DNS data (Rogers 30 ). This model, like the SL model, fully satisfies realizability 

conditions. 


( 2 ) 


Models for the return-to-isotropy term 11^ 
Rotta: 23 


n( 2) = —eCb 


tj 


( 2 . 2 . 6 ) 


where, C = 3.0. T+ 

This model is linear in the Reynolds stress, and contains one model constant, it 

was widely used and adopted in the LRR model. We notice that this model will not 
allow the turbulence to_reach a two-component state, because when any turbulent 
component reduces to q 1 / 9, the model Eq.(2.2.6) will force it to grow. 

Lumlev: 24 

n ( { 2) = -e\/3bij + 7 {b\j + 211^/3)] (2.2.7) 

where, 7 = 0 and 


ft = 2 + — exp(— 7.77 /y/ Re) {72/ y/ He + 80.1 ln[l + 62.4( II + 2.37//)]} 
9 



This model is quasi-linear in the Reynolds stress, because 7 is set to zero, and f3 
is a function of the invariants of Reynolds stress tensor. This model is simple, and 

satisfies realizability. 

Sarkar and Speziale (SS): 25 
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where C\ — 3.4. 

This is a quadratic model in the Reynolds-stress tensor. It satisfies what is 
called the weak realizability condition. Like the Rotta model Eq.(2.2.8), this model 
will not produce unphysical results. However, it will not allow the turbulence to 
approach a two-component state, which could occur in some situations, for example, 
in near-wall turbulence. 

Haworth and Pope (HP ): 26 

ng> = + b% - bl.ib.i + «y/3)]} (2.2.9) 


where C\ = 8.3, Ci = 14.8. 

Eq.(2.2.9) is a slow part of the Haworth and Pope’s model for the situations 
with no mean velocity gradient. This model, like the SS model, will not produce 
unphysical results; however, it will also not allow the turbulence to approach a 
two-component state. 

Choi and Lumlev (CL): 27 
If III > 0, 

nSf = -etybij + 7(6?,- + 2II6ij/3)] (2.2.10.1) 

where, 


(3 = 2 + 


► jpl/2 


7 = 


l + G X 2 

p*F l ! 2 G 


1 + Gx 2 Z 

S = (JJJ/2) 1 / 3 , j] = (-II/3) 1 / 2 

X = -, C = -x 4 + 0.8x 6 
V 

p' = exp [-9.29 / - [296 - 16.2(* + l) 4 ]//} 


Re = , 11= —bijbij/2, III = bijbjkhi/S 

i/C ” 


If m < o, 


ng’ = £9(28) 


(2.2.10.2) 


The model coefficients in Eq.(2.2.10.1) were obtained using their comprehen- 
sive measurements of turbulence relaxing from axisymmetric expansion. Both 
Eq.(2. 2.10.1) and Eq.(2.2.10.2) satisfy realizability; however, Eq.(2.2.10.1) is valid 
only for III > 0, because £ is not defined when III < 0. 

Craft and Launder (CfeL ): 28 


.ng* = -C^bij + 4 c; (bfj - bl^/3)) - lebij ( 2 . 2 . 11 ) 
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where, 


Ci = 3 . 1 (^ 2 > 1 ) 1/2 , C\ = 12 

9 

A<i — 4 bijbji, A3 = 8 bijbjkbki> A 1 — g(-^2 A3) 


This model is tensorially quadratic in the Reynolds stress, and satisfies realizability. 
Yamamoto and Arakawa (YA ): 29 


= -e[ai bij + OL2{b]j - b 2 kk Sij/Z)] 


( 2 . 2 . 12 ) 


Wh6re ’ Q!=2+pF[q ( b 2 kk ) r + \bl k \ a sign{bl k )] 

ctz = 3 (ai — 2) 

p = —12, q = —0.65, r — 0.4, s = 0.45 
F = 1 - ^ b\ k +9 b\ k 

The YA model tried to fit situations with both positive and negative b\ k . It also 
meets the requirement of realizability. 

Shih and Manso ur (SfeM ): 22 

n (2) = _ e {(2.0 + CfF^bn + iKj + (1/3 + 217)&ij + (2.2.13) 


where, 

C f = (1/9) exp(-7.77//Re){72/\/Re + 80.1 ln[l + 62.4(-// + 2.3 III)]} 

T 

1 = 7o(l - F v ), Re - — 

F = 1 +97J + 3//J 
II = — ^b ijbij , III = -bijbjkbki 
£ = 17/20, 7/ = 1/20, 70 = -2 

This model matches the data of Comte-Bellot and Corrsin 21 and meets the require- 
ment that there will be no return to isotropy in the zero Reynolds number limit. 
This model also satisfies realizability. 


Concluding remarks 

We notice that the Reynolds number in all these simulations is low and therefore 
may not represent real turbulence in nature. However, the model terms concerned 
here are mainly pressure related correlations. The fluctuating pressure is no 1 . 
rectly related to the viscosity, hence the pressure related correlation terms may 
not be directly affected by the Reynolds number, especially the rapid term. The 
return-to-isotropy term, which includes the deviatoric part of t e issipa ion ra e 
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tensor may have some dependence on the Reynolds number'. According to the above 
consideration, we think that direct comparisons with low-Reynolds DNS data are le- 
gitimate, although we should keep in mind the possible low-Reynolds number effect 
' of the DNS data. 

We have directly compared five rapid models with fifteen DNS flows: four of 
Rogers et al.’s 30 shear flows, eleven of Lee et al.’s 31 irrotational strain flows (axisym- 
metric contraction, axisymmetric expansion and plane strain). Comparing the per- 
formance of the LRR and SSG models, which are tensorially linear in the Reynolds 
stress, we conclude that the SSG model gives very little improvement over the LRR 
model. In fact in many cases, it is worse than the LRR model. The reason is not 
very clear. However, we notice that the SSG model does not satisfy the normal- 
ization condition which may be a cause for its poor behavior. If we impose this 
constraint on the SSG model, then it will exactly reduce to the LRR model. In fact 
it can be shown that the most general form of the rapid model, which is tensorially 
linear in the Reynolds stress, is the LRR model. Therefore, in general, the treat- 
ment used in the SSG model would hardly give any improvement over the LRR 
model. A natural way to improve the model is to use a more general nonlinear form 
and more general model constraints. A typical example is the SL 21 model. It starts 
with the most general form, using full realizability constraints together with the 
other conventional constraints 34 . The result is a well behaved model. Indeed, from 
the direct comparisons with the DNS data, the SL 21 model and its variation form 
of SM 22 model give the best performance in most of the cases. As to the FLT 20 
model, it is also a nonlinear model. It starts with a tensorially cubic dependance on 
the Reynolds stress with constant coefficients (in general, these coefficients should 
not be restricted to constants). In addition, the two-component conditions of tur- 
bulence have been imposed. However, the FLT model ignores Schwarz’ inequality. 
Its final form contains two undetermined model constants, but one of them is set 
to zero. The performance of the FLT model, from the direct comparisons with the 
DNS data, is better overall than the linear models, but does not compare with the 
performance of the SL and SM models. So from these direct comparisons of the 
rapid models, we conclude that the SL 21 model and its variation form SM 22 are 
clearly the best. Having said this we notice that, as Reynolds 33 pointed out, any of 
these rapid models will not show any effect of rotation on the invariants (II, III ) of 
the anisotropy tensor bij. This is clearly a theoretical deficiency of the current rapid 
models. A further investigation is needed to find out how serious this deficiency 
will be in practice. 

We have directly compared eight return-to-isotropy models with thirty four DNS 
flows: four shear flows and thirty relaxation flows from axisymmetric contraction, 
axisymmetric expansion and plane strain. As was discussed earlier, all the return- 
to-isotropy models are variations of Eq.(2.2.7) derived by Lumley 24 . Therefore 
the differences in the models are due to the different choices of the model coeffi- 
cients. Two linear models are due to Rotta 23 and Lumley 24 (which is quasi-iinear 
in b^). Lumley ’s model satisfies realizability, matches the data of Comte-Bellot 
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and Corrsin 32 and the limit of the final period of decaying turbulence. It per- 
forms perfectly when III < 0. It also compares weU with the DNS data m which 
III > 0 Rotta’s model does not compare with the performance of Lumley s model. 
In fact, the nonlinear models of SS, YA, HP and C&L also do not compare with the 
performance of Lumley’s model. Apparently the model coefficients chosen m these 
models are not appropriate. The CL 27 model is designed for flows with III > 0 and 
is based on their experiments on relaxing turbulence. It does work better than Lum- 
ley’s model when III > 0. Finally, the S&M 22 model is a nonlinear model; it works 
just like Lumley’s model when III < 0. When III > 0, it shows an improvement 
over Lumley’s model according to the DNS data. So from these direct comparisons 
of the return-to-isotropy models, we conclude that the combination of Lumley s 
model and Choi’s model, that is the CL 27 model, will give the best performance. 
The S&M 22 model seems as good as the CL model according to these comparisons. 
Having said this, we notice that the existing return-to-isotropy models do not follow 
the relaxation flows from expansion and plane strain very well. Therefore there is 
still a need to further investigate and improve the returmto-isotropy models. 

For detailed comparison in each flow see the reference 34 . 


2.3 Near-wall behavior of turbulence 

The near-wall behavior of turbulence is re-examined in a way different from that 
proposed by Hanjalic and Launder 35 and followers 36 ’ 37 ’ 38 * 3 - It is shown that at a 
certain distance from the wall, all energetic large eddies will reduce to Kolmogorov 
eddies (the smallest eddies in turbulence). All the important wall parameters, sue 
as friction velocity, viscous length scale, and mean strain rate at the wall, are 
characterized by Kolmogorov microscales. According to this Kolmogorov behavior 
of near-wall turbulence, the turbulence quantities, such as turbulent kinetic energy, 
dissipation rate, etc. at the location where the large eddies become “ K °J™ ogoro 
eddies, can be estimated by using both direct numerical simulation (DNS) data 
and asymptotic analysis of near-wall turbulence. This information will provide 
useful boundary conditions for the turbulent transport equations. As an example, 
the concept is incorporated in the standard k-e model which is then applied to 
channel and boundary layer flows. Using appropriate boundary conditions (based 
on Kolmogorov behavior of near-wall turbulence), there is no need for any wall- 
modification to the k-e equations (including model constants). Results compare 

very well with the DNS and experimental data. xt40A _ U i __ 

Here we only list the results from this study, for the detail see NASA TM 1056 . 

Model equation and boundary condition 

The K-e equations for incompressible flows can be in general modeled as: 

“tt" = \i v + + vrUijWu + Uj t i) — e 

Dt 

^ = [(y + + ClhvrUijWij + 

e 2 

- + WTUi,jkUi,jk 


(2.3.1) 

(2.3.2) 
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Ci = 1.44, C 2 = 1.92 

/i = l, h = 1-0.22 exp (-^), ^ = “7 ( 2 ‘ 3 * 3 ) 

= 1, cr e = 1.3 


These equations are used only for the flow field outside of the turbulence limit point 
where K v is non-zero. Therefore, Eq.(2.3.2) will not have singularity problems 
and will not need any near-wall modifications like other K-e models do. 2,8 
Eddy viscosity: 

v T = C^f^j (2.3.4) 

where, 


C M = 0.09 

fn = [1 - exp(aiR k + a 3 i?£ + a 5 R|)]* 
oi = —1.5 * 10 -4 fl 3 = -1.0*10" 9 a 5 


Rk 


K^y 

V 


-5.0 * 10" 10 


Boundary conditions: at y n = 6i '/u r , 


(2.3.5) 


e v = 0.251^ (2.3.6) 

K v = 0.250u 2 (2.3.7) 


In practical applications, R er and R eoo are large numbers, hence y n /L ( L is the 
length scale of a flow field) is usually very small. Therefore, as an approximation 
we may let y^/L = 0, but e ^ and K v must be given by Eqs.(2.3.6) and (2.3.7) 
respectively. These equations have been applied to the calculations of channel and 
boundary layer flows. 

Comparison of models 

To compare the present model with the DNS data and other models (e.g. Jones 
and Launder 8 , and Chien 2 ), we have made calculations on two channel flows 15,40 
and two boundary layer flows 16 ’ 17 . In the present model, all the model constants 
are the same as used in the standard K-e model 39 . Therefore the present model 
will also be suitable for flows away from the wall. The other two models used here 
for comparison do not have this property. Results are shown in figures 1-4. In 
figure 1 and figure 2, three models are compared with two DNS data for channel 
flows: one with R er — 180, the other with R er = 395. The profiles of mean velocity, 
Reynolds stress, turbulent kinetic energy and its dissipation rate are plotted in these 
figures. The present model is significantly better than the other two models. Figure 
3 shows the similar comparison for a turbulent boundary with R e & = 1410. The 
agreement between the present model and DNS data is excellent. Figure 4 shows 
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the results compared with Klebanoff 17 and other boundary layer experiments. The 
skin friction from DNS data 4 * * * * * * * * * * * 16 is also shown in this figure. The results of present 
model are more consistent with the DNS data than the experiments.^ 

It is also worthwhile to emphasize that the present model equations with the 
standard model coefficients have the simplest form among all two-equation models. 
Hence, we expect that they will have less numerical stiffness in complex turbulent 

flows. 

2.4 Modeling of transition (Z. Yang and T.-H. Shih) 

A model of intermittancy based on the shape factor is added to a two-equation 
k-e model for prediction of boundary layer transition with a free stream turbulence. 
The detailed model and calculations are give by Yang in this year’s Research Briefs. 

2.5 Modeling of compressible turbulence (W.W. Liou and T.H. Shih) 

A two-scale model is proposed based on Hanjalic-Launder’s multiple-scale concept 
for compressible turbulence, in which a distinct scale created by the compressibility 
is modeled separately by considering the effects of pressure-dilatation and dissipa- 
tion dilatation on large-scale energy transfer rate. The detailed model is given by 
Liou in this year’s Research Briefs. 

2.6 Direct numerical simulation of compressible flows (A. Hsu 
and T.-H. Shih) 

In order to have a better understanding of the effect of compressibility on tur- 
bulence, especially the effect of the formation of eddy-shocklets on turbulence, a 
direct numerical simulation of compressible homogeneous shear turbulent flows is 
been performing. The data of all turbulence statistics are very useful for turbulence 
modeling. The detailed simulation is given by Hsu in this year’s Research Briefs. 

3. Future Plans 

Development of second order closure models: pay special attention to the effects 
of inhomogeneity, non-local property, frame-rotation, compressibility, near-wall be- 
havior in the buffer and log-layers. 
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Figure 1. Comparison of models with DNS of 2-D channel with R* r = 
180. 
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Figure 2. Comparison of models with DNS of 2-D channel with R er 
395. 
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Figure 3. Comparison of models with DNS of 2-D boundary layer flow 
with R e e — 1410. 
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